Dealing with many zeros
Count data are often modeled with Poisson or Negative Binomial models.
These assume zeros occur naturally from the count process.
But in real data, we often see too many zeros → zero inflation.
Examples:
Both handle excess zeros but use different assumptions about how zeros are generated.
Zero-Inflated Models
Hurdle Models
Both handle excess zeros but use different assumptions about how zeros are generated.
Zero-Inflated Models
Hurdle Models
Both handle excess zeros but use different assumptions about how zeros are generated:
Zero-Inflated Models
Hurdle Models
Assume two processes:
A binary process → decides if the observation is a certain zero.
A count process → generates counts, including zeros.
\[ \begin{aligned} \text{P}(Y=0) & = \pi + (1-\pi)f(0) & \\ \text{P}(Y=y) & = (1-\pi)f(y)& y=1,2,\dots \end{aligned} \]
Note In zero-inflated models \(f(y)\) is typically discrete (ex Poisson, Binomial, …)
Also have two parts, but different logic:
A binary model determines if the observation crosses the hurdle (i.e. zero vs. positive).
A truncated count model models only positive counts.
Zeros come only from the first process.
\[ \begin{aligned} \text{P}(Y=0) & = \pi & \\ \text{P}(Y=y) & = \frac{f(y)}{1-f(0)}& y=1,2,\dots \end{aligned} \]
Note In hurdle models \(f_Y(y)\) can also be continuous (f.ex Gamma). In these cases: \[ \begin{aligned} \text{P}(Y=0) & = \pi & \\ f_Y(y) & = f(y)& y>0 \end{aligned} \] where \(f(y)\) is scaled so that it integrates to \(1-\pi\)
| Feature | Zero-Inflated | Hurdle |
|---|---|---|
| Zeros come from | Two sources (inflation + count) | One source (hurdle only) |
| Positive part | Includes zeros | Truncated at zero |
| Can use continuous distributions? | ❌ Typically count only | ✅ Yes (Gamma, Lognormal) |
inlabru (Type 1)zeroinflatedpoisson1)zeroinflatedbinomial1)zeroinflatednbinomial1)zeroinflatedbinomial1)To get details about the distributions you can type
The probability of the inflation process \(\pi\) is a hyperparameter therefore cannot be modeled with covariates
nofish persons child count
1 1 1 0 0
2 0 1 0 0
3 0 1 0 0
\[ \begin{eqnarray} \text{Model 1: }\ & Y\sim \text{Poisson}(\lambda)\\ \text{Model 2: }\ & Y\sim \text{NegBinomial}(n,p)\\ \text{Linear predictor: }\ & \eta = \beta_0 + \beta_1x_1 + \beta_2x_2 \end{eqnarray} \]
\[ \begin{eqnarray} \text{Model 1: }\ & Y\sim \text{Poisson}(\lambda)\\ \text{Model 2: }\ & Y\sim \text{NegBinomial}(n,p)\\ \text{Linear predictor: }\ & \eta = \beta_0 + \beta_1x_1 + \beta_2x_2 \end{eqnarray} \]
Negative Binomial distribution
\[ \text{Prob}(Y=y) = \frac{\Gamma(y+n)}{\Gamma(n)\Gamma(y+1)}p^n(1-p)^y \] with \[ E(Y) = \mu = n\frac{1-p}{p} \qquad \text{Var}(Y) = \mu(1+\frac{\mu}{n}) \] and linear predictor \(\eta = log(\mu)\)
The Model \[ \begin{aligned} \text{Model 1 }: y_i|\eta_i & \sim \pi1_{(y=0)} + (1-\pi)\text{Pois}(\lambda(\eta_i))\\ \text{Model 2 }: y_i|\eta_i & \sim \pi1_{(y=0)} + (1-\pi)\text{NegBin}(n(\eta_i),p(\eta_i))\\ \text{Linear Predictor }: \eta_i & = \color{red}{\boxed{\beta_0 + \beta_1 x_{1i} + \beta_2 x_{2i}}} \end{aligned} \]
The code
# define model components
cmp = ~ Intercept(1) + fishing(nofish, model = "linear") +
persons(persons, model = "linear") +
child(child, model = "linear")
# define model predictor
formula = count ~ .
# build the observation model
# Poisson model
lik_pois = bru_obs(formula,
family = "zeroinflatedpoisson1",
data = zinb)
# Negative binomial model
lik_nbin = bru_obs(formula,
family = "zeroinflatednbinomial1",
data = zinb)
# fit the model
fit_pois = bru(cmp, lik_pois)
fit_nbin = bru(cmp, lik_nbin)The Model \[ \begin{aligned} \text{Model 1 }: \color{red}{\boxed{y_i|\eta_i}} & \color{red}{\boxed{\sim \pi1_{(y=0)} + (1-\pi)\text{Pois}(\lambda(\eta_i))}}\\ \text{Model 2 }: y_i|\eta_i & \sim \pi1_{(y=0)} + (1-\pi)\text{NegBin}(n(\eta_i),p(\eta_i))\\ \text{Linear Predictor }: \eta_i & = \beta_0 + \beta_1 x_{1i} + \beta_2 x_{2i} \end{aligned} \]
The code
# define model components
cmp = ~ Intercept(1) + fishing(nofish, model = "linear") +
persons(persons, model = "linear") +
child(child, model = "linear")
# define model predictor
formula = count ~ .
# build the observation model
# Poisson model
lik_pois = bru_obs(formula,
family = "zeroinflatedpoisson1",
data = zinb)
# Negative binomial model
lik_nbin = bru_obs(formula,
family = "zeroinflatednbinomial1",
data = zinb)
# fit the model
fit_pois = bru(cmp, lik_pois)
fit_nbin = bru(cmp, lik_nbin)The Model \[ \begin{aligned} \text{Model 1 }: y_i|\eta_i & \sim \pi1_{(y=0)} + (1-\pi)\text{Pois}(\lambda(\eta_i))\\ \text{Model 2 }: \color{red}{\boxed{y_i|\eta_i}} & \color{red}{\boxed{\sim \pi1_{(y=0)} + (1-\pi)\text{NegBin}(n(\eta_i),p(\eta_i))}}\\ \text{Linear Predictor }: \eta_i & = \beta_0 + \beta_1 x_{1i} + \beta_2 x_{2i} \end{aligned} \]
The code
bru_options_set(control.compute = list(dic = T, waic = T, mlik = T))
# define model components
cmp = ~ Intercept(1) + fishing(nofish, model = "linear") +
persons(persons, model = "linear") +
child(child, model = "linear")
# define model predictor
formula = count ~ .
# build the observation model
# Poisson model
lik_pois = bru_obs(formula,
family = "zeroinflatedpoisson1",
data = zinb)
# Negative binomial model
lik_nbin = bru_obs(formula,
family = "zeroinflatednbinomial1",
data = zinb)
# fit the model
fit_pois = bru(cmp, lik_pois)
fit_nbin = bru(cmp, lik_nbin)[1] "Poisson Model"
mean sd 0.025quant 0.975quant
Intercept 0.22 0.15 -0.08 0.51
fishing -0.90 0.12 -1.13 -0.67
persons 0.65 0.05 0.56 0.74
child -0.95 0.10 -1.15 -0.76
[1] "Negative Binomial Model"
mean sd 0.025quant 0.975quant
Intercept -0.80 0.31 -1.42 -0.18
fishing -0.59 0.25 -1.08 -0.10
persons 0.94 0.12 0.71 1.17
child -1.65 0.19 -2.03 -1.27
Hyperparameters
[1] "Poisson Model"
mean sd
zero-probability parameter for zero-inflated poisson_1 0.48 0.04
[1] "Negative Binomial Model"
mean sd
size for nbinomial_1 zero-inflated observations 0.56 0.10
zero-probability parameter for zero-inflated nbinomial_1 0.06 0.06
Scores
Score Poisson Nbin
1 DIC 1111.36 795.44
2 Mlik -579.97 -424.72
3 WAIC 1207.12 798.86
\[ \begin{eqnarray} \text{Model 1: }\ & P(Y=0) &= \pi + (1-\pi)\ \text{Poisson}(\lambda)\\ \text{Model 2: }\ & P(Y=0)&= \pi + (1-\pi)\ \text{NegBinom}(n,p)\\ \end{eqnarray} \] We check the distribution for one person, non-fisher and no children
df_pred = data.frame(nofish = 0, persons = 1, child = 0)
prob_pois = predict(fit_pois, df_pred,
formula = ~ {
pi0 <- zero_probability_parameter_for_zero_inflated_poisson_1
lambda = exp(Intercept + fishing + persons + child)
yy = 0:20
list(
prob0 = pi0 * c(1,rep(0,20)) + (1-pi0) * dpois(yy, lambda = lambda))
} )
prob_nbin = predict(fit_nbin, df_pred,
formula = ~ {
pi0 <- zero_probability_parameter_for_zero_inflated_nbinomial_1
size = exp(size_for_nbinomial_1_zero_inflated_observations)
mu = exp(Intercept + fishing + persons + child)
p = size/(size + mu)
yy = 0:20
list(prob0 = pi0 * c(1,rep(0,20)) + (1-pi0) * dnbinom(yy,size = size, mu = mu))
})Note: To get the names to input in predict use the bru_names() function:
inlabruThere are two ways to implement hurdle models in inlabru
Strategy 1 Use a modified version of the same distributions defined for the zero-inflated models
Strategy 2 Use a “two-likelihood” trick
Available models (Type 0)
Poisson (zeroinflatedpoisson0)
Binomial (zeroinflatedbinomial0)
Negative Binomial (zeroinflatednbinomial0)
BetaBinomial (zeroinflatedbinomial0)
Advantages
Disadvantages
The Model \[ \begin{aligned} \text{Model 1 }: y_i|\eta_i & \sim \pi1_{(y=0)} + (1-\pi)\text{TruncPois}(\lambda(\eta_i))\\ \text{Model 2 }: y_i|\eta_i & \sim \pi1_{(y=0)} + (1-\pi)\text{TruncNegBin}(n(\eta_i),p(\eta_i))\\ \text{Linear Predictor }: \eta_i & = \beta_0 + \beta_1 x_{1i} + \beta_2 x_{2i} \end{aligned} \]
The code
# define model components
cmp = ~ Intercept(1) + fishing(nofish, model = "linear") +
persons(persons, model = "linear") +
child(child, model = "linear")
# define model predictor
formula = count ~ .
# build the observation model
# Poisson model
lik_pois = bru_obs(formula,
family = "zeroinflatedpoisson0",
data = zinb)
# Negative binomial model
lik_nbin = bru_obs(formula,
family = "zeroinflatednbinomial0",
data = zinb)
# fit the model
fit_pois0 = bru(cmp, lik_pois)
fit_nbin0 = bru(cmp, lik_nbin)Advantages
Disadvantage
The Hurdle Model handles this by splitting the data-generating process into two parts:
This idea can be used in inlabru by using two likelihood
Daily precipitation in Parana state
\[ z_{st} = \left\{ \begin{eqnarray} 0 & \text{ if } &\text{ no rain on day } t \text{ at station } s\\ 1 & \text{ if } &\text{ rain on day } t \text{ at station } s \end{eqnarray} \right. \]
\[ y_{st} = \left\{ \begin{eqnarray} \text{NA} \text{ if } && \text{ no rain on day } t \text{ at station } s\\ \text{rain amout} \text{ if } && \text{ rain on day } t \text{ at station } s \end{eqnarray} \right. \]
\[ z_{st}\sim \text{Binom}(p_{st}, n_{st} = 1); \qquad y_{st}|z_{st}=1\sim\text{Gamma}(a,b) \ \text{with}\ E(y_{st}) = \mu_{st} = a/b \] Where \[ \begin{eqnarray} \eta_t^1 & = \text{logit}(p_t) & = \beta_0^B+u^B(s) \\ \eta_t^2 & = \log(\mu_t)& = \beta_0^G+u^G(s) \end{eqnarray} \]
Implementation
cmp = ~ -1 + Intercept_z(1) + Intercept_y(1) +
space_z(geometry, model = spde) +
space_y(geometry, model = spde) +
local_z(station, model = "iid")
lik_z = bru_obs(formula = z ~ Intercept_z + space_z + local_z,
data = df %>% filter(time==1),
family = "binomial")
lik_y = bru_obs(formula = y ~ Intercept_y + space_y,
data = df %>% filter(time==1),
family = "gamma")
fit = bru(cmp, lik_z, lik_y)Results
We can:
Let’s look again at the fishes examples. We now want to include covariates also in the model for the zero probability:
\[
\begin{eqnarray}
P(Y = 0) = &\pi\\
P(Y=y) = & (1-\pi) \frac{f_Y(y)}{f_Y(0)},& \qquad y =1,2,\dots
\end{eqnarray}
\] We use the same “trick” as before and we now have \[
\begin{eqnarray}
z\sim\text{Binomial}(\pi,n); &\qquad \eta_1 = \text{logit}(\pi) = \beta^1X\\
z\sim\text{Truncated Poisson}(\lambda); &\qquad \eta_2 =\log(\lambda) = \beta^2X\\
\end{eqnarray}
\] In inlabru the truncated Poisson distribution is called nzpoisson
Implementation
cmp = ~ -1 + Intercept_z(1) + Intercept_y(1) +
nofish_z(nofish, model = "linear") + persons_z(persons, model = "linear") +
nofish_y(nofish, model = "linear") + persons_y(persons, model = "linear")
lik_z = bru_obs(formula = z~ Intercept_z + nofish_z + persons_z,
data = zinb,
family = "binomial")
lik_y = bru_obs(formula = y~ Intercept_y + nofish_y + persons_y,
data = zinb,
family = "nzpoisson")
fit = bru(cmp, lik_z, lik_y)Results
mean sd 0.025quant 0.5quant 0.975quant mode kld
Intercept_z 0.70 0.33 0.04 0.70 1.35 0.70 0
nofish_z 0.23 0.28 -0.33 0.23 0.78 0.23 0
persons_z -0.18 0.12 -0.41 -0.18 0.04 -0.18 0
Intercept_y 0.35 0.15 0.05 0.35 0.65 0.35 0
nofish_y -0.83 0.12 -1.07 -0.83 -0.59 -0.83 0
persons_y 0.53 0.04 0.44 0.53 0.61 0.53 0
In inlabru there is no implementation for the truncated binomial (or negative binomial) distribution
We can “trick” inlabru by using one of the implemented Type 0 models \[
y_i|\eta_i \sim \pi1_{(y=0)} + (1-\pi)\frac{f(y)}{f(0)}
\] with a fixed and small \(\pi\) thus creating a truncated distribution!
Implementation
# Define response for binary and truncated neg. binomial
zinb = zinb %>% mutate(z = ifelse(count==0,1,0),
y = ifelse(count==0,NA, count))
# define model components
cmp = ~ -1 + Intercept_z(1) + Intercept_y(1) +
fishing_z(nofish, model = "linear") + fishing_y(nofish, model = "linear") +
persons_z(persons, model = "linear") + persons_y(persons, model = "linear") +
child_z(child, model = "linear") + child_y(child, model = "linear")inlabruThere are mainly two types of models to treat excess of zeros
inlabruThere are mainly two types of models to treat excess of zeros
inlabru (Type1)inlabruThere are mainly two types of models to treat excess of zeros
inlabruinlabruThere are mainly two types of models to treat excess of zeros
inlabru